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We study diffusive magnetic Josephson junctions with four superconducting terminals in the weak proximity 
limit where the leads are arranged in cross form. Employing the linearized Keldysh-Usadel technique, the 
anomalous Green's function and Josephson current are analytically obtained based on a quasiclassical theory 
using the Fourier series method. The derived results may be reduced to non-magnetic junctions by setting the 
exchange field equal to zero. We find that increments of the magnetic barrier thickness may cause a reversal of 
the supercurrent direction flowing into some of the leads, whereas the direction of current-flow remains invariant 
at the others. The reversal direction can be switched by tuning the perpendicular superconducting phases. In 
the non-magnetic case, we find that the supercurrent flowing between the leads in one direction can be tuned 
by changing the superconducting phase difference in the perpendicular direction. These findings suggest the 
possibility of constructing a nano-scale superconducting phase transistor whose core element consists of the 
proposed four-terminal Josephson junction with rich switching aspects. 

PACS numbers: 74.5().+r, 74.45.+C, 74.78.Na 



1. INTRODUCTION 

When a weak link is established between two superconduc- 
tors, a gradient in the superconducting phases can drive a su- 
percurrent through the system. This Josephson effect'"^ and 
the associated current-phase relation in weak links has been 
investigated extensively in previous literature, see for exam- 
ple the comprehensive reviews Refs. 4 and 5 (see also Refs. 6 
and 7 for magnetic Josephson junctions). 

The proximity effect between superconductors and nor- 
mal diffusive metals was first studied by W.L. McMillan in 
1965**. It is known that the electronic properties of a nor- 
mal metal become altered when placed in proximity to a 
host superconductor For instance, the electronic spectrum 
of the normal metal connected to a superconductor exhibits 
a minigap**"'\ Very recently, the key properties of density 
of states (DOS) of a sandwiched normal metal between su- 
perconductors were employed in an experiment for produc- 
ing a superconducting quantum interference proximity tran- 
sistor (SQUIPT)'"^. Moreover, superconductor-normal metal- 
superconductor (S/N/S) Josephson junctions have been stud- 
ied under non-equilibrium conditions where two additional 
normal leads are connected to the sandwiched normal layer It 
has been demonstrated that this type of S/N/S Josephson junc- 
tions is able to produce a tt -junction depending on the applied 
voltage to the normal sandwiched layer'^'"*' . Such tt -junctions 
may also be observed in three terminal junctions'^'"*^. 

The proximity-induced interplay between superconductiv- 
ity and ferromagnetism in hybrid structures is also known to 
establish intriguing physical phenomena. The wavefunction 
describing the leakage of Cooper pairs inside a ferromagnet 
oscillates in a damped fashion. One of the most interesting 
phenomena in the proximity of ferromagnetism and supercon- 
ductivity is O-TT transition which may occur in superconductor- 
ferromagnet-superconductor (S/F/S) junctions*'-''"^^. The 
transition usually occurs over a narrow length ^i? = \J Dp jh 
in which Dp and h represent the diffusion constant and the 
exchange field of the sandwiched ferromagnetic layer, respec- 



tively. At this crossover point, the minimum energy of junc- 
tion is switched between zero and vr-superconducting phase 
difference by changing the energy scales of the system such 
as Thouless energy, exchange field and temperature. Also it 
has been demonstrated that the spin-flip scattering may render 
the junction energy minimum from to tt^-^^-^i ^jjjJ ^.^y j-jjg 
supercurrent itself may become spin-polarized if the magneti- 
zation texture is inhomogeneous-'^. 

So far in the literature, the main emphasis has mostly been 
on one-dimensional systems where two superconductors are 
coupled via e.g. a constriction or diffusive metal. On the 
other hand, the interplay between multiple superconducting 
terminals'^ in a Josephson junction would require an exten- 
sion to higher dimensions^*'"'''. This in turn complicates the 
analytical treatment of the system, and one is usually forced 
to resort to numerical means within the diffusive regime'*". 
It would therefore be of interest to clarify how the transport 
characteristics of a diffusive ferromagnetic Josephson junc- 
tion is influenced by the presence of multiple superconducting 
phase differences, and also to provide an analytical framework 
for studying such phenomena. Multi-terminal Josephson point 
contacts had intensively been investigated (both AC and DC 
characteristics) using the Ginzburg-Landau theory^^'^^'^ and 
was followed by studying the four-terminal S/N/S Josephson 
junctions in the clean limit via the Eilenberger equations'"'"^**. 
Interesting phenomena such as phase dragging (the produc- 
tion of phase difference between two terminals by means of 
phase variation between other terminals), magnetic flux trans- 
fer and bistable states were found due to non-local coupling 
and additional degrees of freedom in such classes of Joseph- 
son junctions^''''''^'''^. Such point contacts also have been fab- 
ricated and intensively studied in experiments^'*. 

Motivated by this, we consider in this paper a diffusive 
Josephson junction with four superconducting leads where are 
arranged in a cruciate form and study the supercurrent flowing 
in this junction. The superconducting leads are separated by a 
metal that may or may not be ferromagnetic. We use the quasi- 
classical Usadel equations in the diffusive regime and formu- 
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late the current-phase relation as a function of all the available 
parameters in the system such as superconducting phases in 
the magnetic junction. We recover the results of Refs. 36 and 
37 obtained in the clean S/N/S junctions: namely, when the di- 
mensions L (length) and W (width) of the sandwiched metal 
are comparable to each other, i.e. L ~ W , the standard sinu- 
soidal supercurrent is strongly modified by all the condensate 
phases. We also use a phenomenological Ginzburg-Landau 
theory to confirm our analytical expressions obtained via the 
quasiclassical framework. In particular, we demonstrate that 
the Josephson current flowing between leads along one axis 
may be tuned via the superconducting phase gradient in the 
perpendicular direction. 

Moreover, we find that increments of the magnetic bar- 
rier thickness may cause a reversal of the supercurrent direc- 
tion flowing into some of the leads, whereas the direction of 
current-flow remains invariant at the others. These findings 
are suggestive in terms of designing a nano-scale supercon- 
ducting phase transistor where current switching effects in one 
direction is possible by variation of macroscopic supercon- 
ducting phase in the perpendicular direction as has also been 
pointed out in Ref 38 and 39 for ballistic contacts. 

The paper is organized as follows. In Sec. 2 we present our 
main analytical findings. In Subsect. 2.1 the basic equations 
of the quasiclassical method are presented and in Subsect. 2.2 
the cruciate Josephson junction is studied analytically via the 
Green's function method. We formulate the current-phase re- 
lation as a function of the four superconducting phases for 
a magnetic Josephson junction. In Subsect. 2.3 we confirm 
our results and findings via a macroscopic Ginzburg-Landau 
theory. In Sec. 4 we employ a 'Jacobi' numerical method"*^ 
(which shall be explained in detail) and investigate the behav- 
ior of the supercurrent which confirms our analytical derived 
expressions in Subsect. 2.2 and their dependencies on the su- 
perconducting U(l) phases, also the behavior of junction is 
analyzed in more detail. Sec. 4 is devoted to the study of 
the supercurrent behavior in S/F/S four- terminal junctions as 
a function of ferromagnetic barrier thickness. Concluding re- 
marks are finally given in Sec. 5. 



2. THEORY AND ANALYTICAL DISCUSSIONS 

We consider four superconducting leads coupled via a fer- 
romagnetic or normal diffusive metal. As in Fig. 1, the nano- 
scale diffusive metal is assumed to be located in the xy plane, 
where x G [0,L] and y £ [0, M^]. The four superconduct- 
ing terminals are assumed to have equal magnitudes for the 
gap A and are connected to each edge of the diffusive strip. 
The suppression of the pair potential is neglected near inter- 
faces due to a low interface transparency and the supercon- 
ducting phases are assumed to be different in each of the four 
terminals: 6'up, 6'downi ^left and 6'right. One may expect that su- 
perconducting correlations inside the system interfere, result- 
ing in a quite complicated coherent system. The S/F/S sys- 
tem is studied in the diffusive limit and current-phase rela- 
tionship is obtained at each terminal similar to clean S/N/S 
four-terminal junctions''*''^^. In our approach, we start with a 




FIG. 1: Experimental schematic setup of the cruciate Josephson 
junction. The junction is assumed to lie in the xy plane with inter- 
faces located at a; = 0, L and y = 0, W. The four spin-singlet su- 
perconductors have different superconducting phases: Oup, Odown, dicft 
and bright- Exchange field h, is assumed to be oriented in the z direc- 
tion perpendicular to the sandwiched layer plane. 

magnetic four-terminal Josephson junction and derive our an- 
alytical results for the magnetic system. We then may achieve 
the non-magnetic Josephson junction characteristics by set- 
ting the magnetic exchange field h equal to zero. 



2.1. Microscopic Green's function approach 

In this subsection, we present basic equations of the qua- 
siclassical Keldysh-Usadel method. In order to study the 
transport properties of the proposed four-terminal device, we 
employ the quasiclassical method. In the diffusive regime, 
due to the existence of strong scattering sources, quasiparti- 
cles' momentums are integrated over all directions in space. 
In this case, the Eilenberger equations reduce to the Usadel 
equations'^. Under equilibrium conditions, the system under 
consideration can be described by a 4 x 4 matrix propagator 
in Nambu space: the retarded Green's function . The total 
Green's function describing the system compactly reads'*^: 

where the meaning of the .7. -operation depends on the nota- 
tion adopted. In our notation, it denotes complex conjuga- 
tion and a change in sign for the energy argument. The ad- 
vanced and Keldysh blocks are made from retarded block by 
G^ = -(TgG^rg)^ andG^ = tanh(/3£) (G^-G^^) in which 
T3 is the Pauli matrix and /? = ksT /2. In the presence of 
exchange energy h — {hx, hy,hz) inside the ferromagnetic 
layer, the Usadel equation can be give by; 

D[d, G[d, G]] + i[ep3 + diag[h • a, (h • a)^, G] = 0, (2) 

where and a; are 4 x 4 and 2x2 Pauli matrixes, respectively. 
Here D is diffusive constant of the sandwiched medium. Also, 
£ is the quasiparticles' energy which is measured from Fermi 
surface. 

The so-called weak proximity regime occurs in the case 
of very low transparent interfaces or for temperatures near to 
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the critical temperature of the superconducting leads. The su- 
perconducting correlations leak into the ferromagnetic region 
weakly and so the normal and anomalous Green's functions 
can be approximated by g ~ 1 and / <C 1, respectively. In 
this limit one can linearize the Usadel equation which yields 
a set of uncoupled complex boundary value partial differen- 
tial equations. The energy representation is used in this paper, 
however, one may reach the Matsubara representation by re- 
placing £ iujn, where a;„ = (2n + l)TrkBT are Matsubara 
frequencies. For the sake of simplicity, a uniform exchange 
field for the ferromagnetic layer is considered throughout the 
paper i.e. h = (0, 0, — h). In the weak proximity regime 
that mentioned above, the Green's function read^' 



(3) 



in fact, we have expanded the Green's function around the 
bulk solution Go as G ~ Go + /, where Gq = diag(l_,-l_)^. 
The retarded Green's function now can be given by; 

/I /f(e)\ 

1 f^ie) 

[-ffi-e)r -1 





(4) 
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If we assume that the exchange field is uniform throughout 
the sample and is oriented in the z direction, so the Usadel 
equations reduce to two dimensional form as belows: 

Olfii-s) + dlf^i^s)-^-^^^fS{^s) = 0, (5) 

dl[f^{e)r+dl[f^{e)r - ^-^^^meW = 0. (6) 

We employ the Kupriyanov-Lukichev boundary conditions at 
F/S interfaces'*^ and control their opacities using a parameter ( 
that depends on the resistance of the interface and the diffusive 
normal region; 



C(G9G)-n = [Gbcs(0),G], 



(7) 



where n is a unit vector denoting the perpendicular direction 
to an interface. The bulk solution , Gbcs for a s-wave super- 
conductor is'**; 



r<Ti. 
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i?(e) = arctanh(J-^), 
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s(e) = sinh(i?(e))e*' 



-A 



sgn(e) 



. - A2 
c(e) = cosh(z9(e)) 
le' 



e(£2 - A2) - 



= e(A2-e2) 



- A- 



:e(£2_ a2) 



VA2 - e- 



A is superconducting gap in the s-wave superconductors and 
the Heaviside step-function is denoted by 8(e). In this paper, 
we have defined 6^^, O^t, 6i, 9^ as the condensate phases in the 
up, down, left and right superconductor leads, respectively. If 
we now open up the compacted boundary conditions Eq. (7) 
at left F/S interface for instance, a; = 0, we reach at; 

(Ca,-c* (£))[/«(£)]* =Ts*(e)e-^^', (9) 



and at x 



(Ca, + c* (£))/«(-£) = Ts*(e)e^ 
(Ca, + c*(e))[/^(£)]*=±s*(e)e 



(10) 



Also aty — 

{<:dy-c*{emi{e)r ^Ts*{e)e-^'^, (11) 
and aty ~ W the boundary condition takes the below form 

{Cdy + c*{e))f^i-e) = Ts*ie)e''- 

iCdy + c*{e))[fS{e)r ^±s*{s)e-^'-. (12) 

In the equilibrium conditions, the current density vector is 
given by Keldysh block as 



J(R) = Jo / deTr{p3(G[9,G])^} 



(13) 



here Jo is a normalization constant. The current density vector 
determines the direction and amplitude of current density in- 
side the sandwiched layer as a function of coordinates. If we 
substitute the total Green's function Eq. (1) into the current 
density relation namely, Eq. (13) we arrive at: 

J(R) = Jo r detanh(£/3) {/i^(~£)V[/f (e)]* 

+/^(-£)V[/«(e)]* - f^isMf^i-eW - f^ie) 
V[/f (-£)]* + [f"i-s)]*^f^ie) + [/f (-e)]*V/i^(e) 
-[/f (£)]*V/i^(-e) - [/i^(e)]*V/f (-e)} . (14) 

To obtain total supercurrent flowing through the junction, for 
example at right superconducting gate, one needs to perform 
an integration of Eq. (13) over the y coordinate , 7(0) — 

loj !dydeTr{pMd,g\)''}- 

At this point it suffices that Eqs. (5) be solved together with 
appropriate boundary conditions (i.e. Eqs. (9), (10), (11) and 
(12)) in order to capture the transport characteristics of the 
present class of Josephson junctions in the diffusive limit. 

2.2. Analytical microscopic discussions 

In this subsection we derive explicit analytical expressions 
describing the supercurrent at each superconducting terminal. 
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To this end, we consider the weak proximity limit of diffu- 
sive regime where the Keldysh-Usadel method yields a set of 
uncoupled complex elliptic partial differential equations. The 
simplified Usadel equations and corresponding boundary con- 
ditions are give by Eqs. (5), (6), (9), (10), (11) and (9). For 
simplicity in our analytical calculations we exclude first-order 
terms of the anomalous Green's function in the Kupryianov- 



Lukichev boundary conditions, Eq. (7). We use the Fourier 
series method in the presence of non-homogenous boundary 
conditions and obtain analytical solutions for the Usadel equa- 
tions. The method leads a somewhat lengthy solutions, for in- 
stance one of the anomalous components of Green's function 
namely, /+ (e) after long calculations is given by Eq. (15); 
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W(£ + fe)(-l)'cOs(j^) . 



(15) 



The length and width of the ferromagnetic region sandwiched 
between the superconductors are denoted by L and W. As can 
be seen, the anomalous component of the retarded Green's 
function depends on all four condensation phases, which in 
turn leads to an interference between these superconducting 
phases in the Josephson current. In Eq. (14) there are 8 dif- 
ferent terms of anomalous component of Green's function in- 
volved the supercurrent relation. Therefore, one must find 8 
similar solutions as Eq. (15) for other terms and substitute 
them into the supercurrent relation Eq. (14) in order to ob- 
tain the supercurrent at one terminal. To obtain analytical 
solutions for the total supercurrent flowing at the other su- 
perconducting terminals, one must repeat the latter described 
process. We have done so and arrived at the analytical ex- 
pressions describing the supercurrent in the system as follows. 
Supercurrent at x = 0,L terminals are obtained as 



4(j; = 0) _ f°° de A2tanh(/3e) ^ T WD 



mo 



LWCie + ah) LWC^{e + ah) 



) sm{9i — Or 



and also at the W = 0,L terminals: 

lyiy = 0) d£ A2tanh(/3£) 



(16) 



lAx L) 
In 



de A2 tanh(/3e) 

a;; A2 - £2 



E ( 



WD 



■L^e{e + ah) 



SWD^ (~l)^(£ + afe) 
+ W h DW/L^+4{e + ah)^>''''^^' ~ ^'■>' 

D sm{9d - 9r) D sin(6'„ - dr) 
LWQ^{e + ah) ^ LWC,^{e + ah) 



(17) 



_^Ao A2-£2 



cr=± 



8LD^ i-iy{e + ah) 

W^^C f-^ D^IH'^ /W^ + 4(£ + ah)^ 

Dsmi9d-er) Dsinjed-ei) 
LWC (e + (^h) LWC^ (e + ah) 



LD 



■W^C{e + ah) 



(18) 



lyiy = w) 



lo 
8LD 



de A2 tanh(/3£) 



E 



-coAo A2 



i-lYie + ah) 



El 



LD 



W^C f-^ D'HH'i/W^ + 4(£ + cr/i)2 • 

D sm{ei -6u) D sm{0r - Ou) 
LWe{e + ah) ^ LWC^ie + ah) 



W^C{e + ah) 
siniOd -0u) + 

(19) 



a — ± comes from the spin-dependent nature of the ferro- 
magnetic material which is sandwiched between the four su- 
perconducting terminals. To be more specific, Ix{x — 0), 
Ix{x ~ L), Iy{y — 0) and Iy{y — W) represent the Joseph- 
son current in the x direction at x — 0,L and y direction at 
y = 0, W, respectively. The above currents involve three si- 
nusoidal terms whose arguments include phase differences of 
the lead which supercurrent is being calculated at and the three 
other terminals. As expected, the obtained supercurrents show 
explicitly that this interfering terms in the x and y directions 
vanish for large L and W, respectively. This fact is also found 
in ballistic junctions^*'^^^. In these two limits, either large L 
or W, the system takes on quasi-one dimensional features and 
we recover the well-known standard sinusoidal Josephson re- 
lation for the supercurrent. However, in the opposite regime 
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where L ~ W, the proximity-induced order parameters from 
the superconducting terminals overlap substantially and ad- 
ditional terms compared to the one dimensional case appear 
in the expressions for the supercurrent. As we shall see, the 
supercurrent can behave strongly different from one dimen- 
sional junctions as a function of the phase in one supercon- 
ducting terminal due to this overlap. In fact, the supercurrent 
is a function of a superposition of sinusoidal phase differences 
between the different superconducting leads and one may ex- 
press the supercurrent relations as I{xi) — J2j Ij sin{di — dj) 
in weakly coupled systems'*^ The conservation of charge 
current is also satisfied by the current relationships namely, 
Eqs. (16), (17), (18) and (19). It can be verified explicitly 
that; 

I,ix = 0) + lyiy = 0) = h{x = L)+ Iy{y = W). (20) 

which constitutes the Kirchhoff law of electricity. We will 
proceed to investigate and justify the obtained analytical su- 
percurrent numerically and study how they depend on the su- 
perconducting phases of the terminals. First, we compare our 
analytical expressions for the supercurrent with the results ob- 
tained via a macroscopic Ginzburg-Landau theory in the next 
subsection. 



2.3. Ginzburg-Landau approach: analytical macroscopic 
discussions 

In this subsection, we make a complementary discussion 
and examine qualitatively the quasiclassical findings of the 
previous subsection by comparison with a phenomenological 
Ginzburg-Landau (GL) theory"". The phenomenological ap- 
proach is a macroscopic theory which is unable to explain the 
microscopic mechanism underlying superconductivity, but in- 
stead describes the macroscopic properties near a phase tran- 
sition of the system by writing the free energy as an expan- 
sion in the order parameter We note that the smallness of the 
superconducting order parameter may be compared directly 
with the weak proximity effect regime in the quasiclassical 
theory for temperatures near Tc- We assume here that the 
normal regions characteristic length scale (d) satisfies ^ ^ d 
where ^ is the coherence length. In this case the condensa- 
tion wavefunctions overlap effectively via the proximity ef- 
fect. It is instructive to briefly consider first the one dimen- 
sional case, where one may write an ansatz for the wavefunc- 
tion as follows^'^': 

^ = yjie'^^X + ij2e''^^l-X). (21) 

Here, ipj is the amplitude of the condensate wavefunction in 
region j — 1,2 while 9j is the corresponding superconduct- 
ing phase. The function X is unknown, but assumed to satisfy 
A" — > 1 inside region 1 while X —i' inside region 2. We 
now generalize this ansatz to the present four-terminal two di- 
mensional case. Assume that deep inside the superconducting 
banks the order parameter is given as 

= V„e*«" , Vde*^^ , Vze*'' , Ae'''' ■ (22) 



Inside the contact region, the four condensation's wavefunc- 
tions overlap and consequently we expect a solution as 

= ^Pre''^xy{i -y) + ^bie^'' (1 - x)yii -y) + 

Ae''^-yX{l -X)+ i^de'^^l - y)X{\ - X), (23) 

here we have generalized the mentioned one dimensional 
ansatz for the four-terminal junction. The functions X and 
y satisfy the following asymptotic behavior: A" — > in the 
left. A" -> 1 in the right, — > in the bottom and 3^ — > 1 
in the top superconductors. The supercurrent density can now 
be defined by the second GL equation"*-^^ 

L = ^Im{V'*VV'}, (24) 

pm 

where a and f3 are phenomenological coefficients in the GL 
theory. After some calculations, we find the following expres- 
sions for ]x and ]y, the supercurrent components in the x and 
y directions, 

ice = - y)y {-y{i - y)A^Pr sin(0, - 9r)~ 

- y)AAsHOd - Or) ~ x^yi^^4,rsm{e^ - e^) 

+ (1 - Xf{l - y)A'iPi smiOd - 0i) + 
y{l-X)^i;uAsm{eu-ei)} (25) 
jy = y (1 - X)X{-X{1 - X)AA sin(0„ - Od)- 

y (1 _ X)i;iil;^ sin(0, - 9^) - y'^X^^ sin(0, - 0„) 

+(1 - yf(\ - x)AA siniOi - ed)+ 
x{i-y)^i'dAM0r-9d)} (26) 

in which the prime sign denotes derivation. The obtained 
results illustrate that, for instance in the terms coupling 
the top and bottom superconducting terminals vanish. In this 
way, we see that the phenomenological GL approach produces 
identical dependencies on the superconducting phase differ- 
ences as the microscopic approach using quasiclassical the- 
ory. Direct comparison with e.g. Eqs. (16) and (17) in the 
appropriate hmits for X shows consistency with Eq. (25). 



3. FOUR TERMINAL NON-MAGNETIC JOSEPHSON 
JUNCTION 

In this section, we first set /i = (the exchange field of 
ferromagnetic layer) and consider an S/N/S junction. Basi- 
cally, there are two methods for inducing a supercurrent into 
our Josephson system: 1) via an external flux where the exter- 
nal magnetic field penetrates the junction through a SQUID- 
like geometry and 2) via a current-bias where the supercur- 
rent is injected into the system. A combination of these two 
methods is also possible by utilizing different configurations 
of a multi-terminal system (for a comprehensive investiga- 
tion of such possibilities, see Refs. 33, 34, 35, 36, 37). The 
supercurrent at each terminal can be generally expressed as 
li = J2i j li.j sin(6'i — 6j). Thus if one is able to tune the su- 
perconducting phases independently, the supercurrent will be 
a 27r-periodic function of one of the superconducting phases. 
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3.1. Numerical justification of current phase relationships 

In this subsection, we discuss the analytical findings ob- 
tained in the previous section and present numerical results 
using a real energy representation. In the actual plots, we 
consider a temperature T = O.OSTc and also set the normal 
region's length and width to i = ~ 2.5^5. In this repre- 
sentation, we normalize lengths against ^5 and introduce the 
Thouless energy Sx = [hD/L^). Also, we have normalized 
the quasiparticles' energy by the superconducting gap at zero 
temperature Aq and consider units so that h~ ks = 1- More- 
over, we add a small imaginary number ry/Ao — 0.1 to the 
quasiparticle energy to account for inelastic scattering which 
leads to a finite lifetime for quasiparticle excitations. Setting 
C, —1 ensures the validity of weak proximity in numerical cal- 
culations. Solving numerically the resultant complex bound- 
ary value partial differential equations, the approximate so- 
lution components of the Usadel equation are assumed to be 
linear combinations of bicubic Hermite basis functions, and 
required to satisfy the Usadel equations (5) and (6) exactly at 4 
collocation points in each subrectangle of a grid, and to satisfy 
the boundary conditions exactly at certain boundary colloca- 
tion points. We mention in passing that we include first-order 
terms of the anomalous Green's function in the Kupryianov- 
Lukichev boundary conditions, as done in Ref. 22, in contrast 
to the usual approximation in the literature where such terms 
are discarded. By doing so, we improve the accuracy of the 
analytical solution in our numerical investigations. Finally, 
the linear algebraic equations resulting from the collocation 
method, which are highly nonsymmetric and thus difficult to 
solve using iterative and sparse direct solvers, are solved using 
a "Jacobi" conjugate-gradient method, which means that the 
conjugate gradient method (Section 4.8 of Ref. 24) is applied 
to the preconditioned equations D^^A^Ax — D^^A^h, 
where D is the diagonal part of A'^A. For a generalized dis- 
cussion see Ref. 23. The same framework was very recently 
used in Ref. 42 to study the anomalous Fraunhofer pattern 
appearing in an inhomogeneous S/F/S structure. 

In order to clarify the behavior of the supercurrent in the 
present four-terminal Josephson junction with respect to con- 
densate phases of the four superconductors, we use the fol- 
lowing strategy. We focus on the behavior of the supercurrent 
with respect to one superconductor's phase (the left one) and 
set two phases equal to zero: Ofjoy^^ = Bright — 0, while vary- 
ing ^up. The motivation for this is to see if the supercurrent 
flowing in one direction can be tuned explicitly by the super- 
conducting phase difference in the transverse direction, which 
would correspond to a superconducting phase transistor-like 
device. 

In general, the supercurrent inside the normal diffusive re- 
gion is described by a vector field and depends on the position. 
The total flowing current is conserved, as we have proven an- 
alytically. We focus here on the supercurrent flowing into and 
out of the terminals, i.e. at the positions x ^ 0, y ^ 0, x = L 
and y ^ W gates. The results are shown in Fig. 2 where we 
plot the supercurrent at the four gates as a function of left su- 
perconducting phase where 0„ is varied while Oct — 6^. — 0. 
The top left frame shows the supercurrent at x = as a func- 




FIG. 2: Top left: Supercurrent in the x direction as a function of left 
condensation phase, 9\at, at left superconductor gate i.e. 2; = 0. Top 
right: Supercurrent in the x direction vs left superconducting phase 
,Sicft, at right superconductor gate i.e. x = L. Bottom left: Super- 
current in the y direction as a function of left condensation phase at 
down superconductor gate i.e. y = 0. Bottom right: Supercurrent in 
the y direction vs left superconducting phase at up gate i.e. y — W. 
Here other superconductor phases namely, 9^^ and 6'down are assumed 
to be zero. 



tion of the left superconducting phase, top right is the super- 
cuiTent at X = L, bottom left frame displays the supercurrent 
at y — 0, and finally the bottom right frame shows the su- 
percurrent at y = W. The standard sinusoidal current-phase 
relation appears at all gates in the special case where 9u is 
equal to zero. This behavior can be understood by consider- 
ing Eqs. (16), (17), (18) and (19). In this case, only terms with 
sin(6';) survive and the supercuiTent exhibits a pure sinusoidal 
relation vs 9i. When 9^ increases, the phase shift effectively 
adds a constant which can be either positive or negative. In 
particular, the cuiTents at a; = L and y — shift either up- 
wards or downwards depending on the value of as can be 
understood by looking at Eqs. (17) and (18): a change in 9^ 
only varies constant terms involving sin(6'„). 

In contrast, variation in 0„ influences the currents at x — 
and y = W in a more complicated manner. In this 
case, there is an explicit dependence on the phase difference 
9i — 9u, which induces a strongly non-sinusoidal behavior in 
the current-phase relation. Interestingly, we see that it is pos- 
sible to cancel out the current even for a finite value of 9i 
by choosing 9^ appropriately. This observation suggests that 
the present four-terminal device can act as a superconducting 
phase transistor where the phase difference in one direction 
controls the supercurrent flowing in the perpendicular direc- 
tion. The underlying mechanism behind this is the interfer- 
ence between the condensate wavefunctions in the diffusive 
normal region, which results in an intricate phase-dependence 
of the supercurrent as shown in the analytical results. 
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FIG. 3: Critical supercurrent as a function of the normalized junc- 
tion length L/^s at different superconducting gates and for various 
values of 6up, the superconducting phase of the up terminal. Top 
left: at the left superconductor gate i.e. a; = 0. Top right: at the right 
superconductor gate i.e. x = L. Bottom left: at the down supercon- 
ductor gate i.e. y = 0. Bottom right: at the up gate i.e. y = W . The 
other superconductor phases are fixed at zero. 



4. FOUR-TERMINAL MAGNETIC JOSEPHSON 
JUNCTION 

In this section, we consider a four-terminal Josephson junc- 
tion with a ferromagnetic barrier where the exchange field of 
the magnetic layer is oriented along the z direction. In the 
usual two-terminal magnetic Josephson junctions, an incre- 
ment of the ferromagnetic barrier thickness not only reverses 
the current direction at particular thicknesses but also renders 
the minimum of junction energy to change from supercon- 
ducting phase difference to a tt phase. The phenomenon is 
so called O-tt transition. As has been discussed in Ref. 36 the 
junction energy where there are several superconducting leads 
can be expressed as Ej — J2j<i 7j,i(l ~ cos{6j — 6i)). Here, 
the i and j indices stand for the ith and jth superconducting 
leads. Below, we demonstrate that an increment in the thick- 
ness of the ferromagnet can reverse the flow of supercurrent 
into a pair of the superconducting terminals (along the direc- 
tion of increment), whereas the current direction in the other 
terminal pair remains unaltered. 



4.1. The behavior of critical supercurrent as a function of 
magnetic barrier thickness 

We here present a numerical study of the transport proper- 
ties of four-terminal ferromagnetic Josephson junctions. Al- 
though the numerical results are confirmed by the analyti- 
cal expressions presented in Sec. 3, we include first-order 
terms of the anomalous Green's function in the Kupryianov- 
Lukichev boundary conditions in contrast to the approxima- 



tion used for deriving the analytical expressions for supercur- 
rent where such terms are dropped. We now consider a non- 
zero value of the ferromagnetic exchange field h. For a weak, 
diffusive ferromagnetic alloy such as Pd^jNii-j;, the exchange 
field H/Aq is tunable by means of the doping level x to take 
values in the range meV to tens of meV. Here, we will fix 
h = IOAq, which typically places the exchange field h in the 
range 10-20 meV. In order to investigate the effects of mag- 
netic barrier thickness on the supercurrent at each terminal 
and the influence of the various superconducting phases, we 
follow a similar strategy as in the previous section. 9i is var- 
ied from to 2tt where magnetic barrier length, L, is being 
varied from L = to L = 5^5. The other superconducting 
phases are fixed at zero except 6'„ which is changed in order 
to demonstrate the possible influence of the other supercon- 
ducting phases. The critical value of the supercurrent at each 
terminal is calculated separately for each value of 

Fig. 3 indicates the behavior of critical supercurrent at 
each superconductor lead as a function of normalized junc- 
tion length for various values of 6'„. The top left frame 
exhibits the critical current at left terminal. Except for 0„ = tt 
which shows two points changing the supercurrent direction, 
the other values give rise to one sign-change in the critical cur- 
rent. Identical qualitative behavior appears for the current at 
the right terminal except when 6^ = 0, as shown in the bottom 
left frame. Top and bottom right frames exhibit the critical su- 
percurrent vs at the down and up terminals, respectively. 
The critical supercurrent at the two terminals show a smooth 
function of which is in stark contrast with the behavior 
of the critical supercurrent at the left and right terminal. Thus, 
the increment of the junction length primarily affects the crit- 
ical supercurrent flowing into leads along the same direction 
of the increment. Moreover, the direction of the current can 
be drastically switched by tuning the superconducting phase 
of up terminal. In contrast, the current flowing into the super- 
conducting banks perpendicular direction to junction length 
increment is left unchanged. This class of multi-terminal fer- 
romagnet Josephson junction then offers an interesting syn- 
thesis between and 7r-states, and possibly ^-states, due to 
the fact that the coefficients Ij can change sign depending on 
the junction parameters such as L and W. 



5. CONCLUSIONS 

In conclusion, we have studied a four-terminal Josephson 
junction where a diffusive normal or ferromagnetic metal with 
sides L and W is sandwiched among four s-wave supercon- 
ductor leads. We have obtained explicit analytical results 
using the quasiclassical Keldysh-Usadel method for the su- 
percurrent in the system. We find that the wavefunctions of 
the four superconductors interfere efficiently when L ^ W 
and modifies the standard sinusoidal current-phase relation 
which confirm previous findings in ballistic junctions. These 
findings are confirmed qualitatively by using a macroscopic 
Ginzburg-Landau theory. We have presented numerical re- 
sults for the behavior of the supercurrent, and demonstrated 
that the current flowing along one axis may be tuned by the 
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superconducting phase-difference along the perpendicular di- 
rection. It is demonstrated that such four-terminal junctions 
can provide a rich switching circuit element (due to additional 
degrees of freedom in comparison with one-dimensional two- 
terminal Josephson junctions) where the various supercon- 
ducting phases influence considerably the current behavior 
at the terminals. In particular, we show that a reversal in 
critical current direction as a function of junction length can 
be strongly switched by means of variation of superconduct- 



ing phase of perpendicular terminals. The present investiga- 
tions of diffusive cruciate magnetic Josephson junction may 
provide new perspectives for the design of a superconduct- 
ing phase switches where can be used in quantum circuits as 
switching elements. 
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